clear all
set mem 4000m
set matsize 2000
set maxvar 8000

use build.dta


********************************
* NEW VARIABLES
********************************

* Per capita wage in manufacturing, in real dollars
* Manufacturing payroll of production workers over average number of production workers 
g wage0    = (pcmwage00)/(cpi0/100)
g wage20   = (pcmwage20)/(cpi20/100)
g wage30   = (pcmwage30)/(cpi30/100)
g wage40   = (pcmwage39)/(cpi40/100)
g wage50   = (pcmwage47)/(cpi50/100)
g wage60   = (pcmwage58)/(cpi60/100)
g wage70   = (pcmwage67)/(cpi70/100)
g wage80   = (pcmwage82)/(cpi80/100)
g wage90   = (pcmwage87)/(cpi90/100)
g wage2000 = (pcmwage97)/(cpi2000/100)


*use lagwage measures to avoid mechanical correlation between L and w
replace wage2000=lagwage00 //92 CBP
replace wage90=lagwage90 //82 EC
replace wage80=lagwage80 //72 EC
replace wage70=lagwage70 //62 EC
replace wage60=lagwage60 //52 EC
replace wage50=sqrt(lagwage60*lagwage30) //interpolate
replace wage40=lagwage30  //this comes from 1930 Decennial



* Per capita wage in trade, in real dollars
* total payroll in wholesale establishments + retail establishments / (total workers in retail+ in wholsale)
g twage30  = pctwage30/(cpi30/100) 
g twage40  = pctwage40/(cpi40/100)
g twage50  = pctwage54/(cpi50/100)
g twage60  = pctwage63/(cpi60/100)
g twage70  = pctwage72/(cpi70/100)
g twage80  = pctwage82/(cpi80/100)
g twage90  = pctwage87/(cpi90/100)
g twage2000  = pctwage97/(cpi2000/100)

* Agricultural values

g lnfaval890 = ln(faval890/(cpi890/100))
g lnfaval0 = ln(faval900/(cpi0/100))
g lnfaval10 = ln(faval910/(cpi10/100))
g lnfaval20 = ln(faval920/(cpi20/100))
g lnfaval30 = ln(faval930/(cpi30/100))
g lnfaval40 = ln(faval940/(cpi40/100))
g lnfaval60 = ln(faval959/(cpi60/100))
g lnfaval2000=ln(faval2002/(cpi2000/100))


* Median family income
gen lnmedfaminc50  = ln(medfaminc50/(cpi50/100)) 
gen lnmedfaminc60  = ln(medfaminc60/(cpi60/100)) 
gen lnmedfaminc2000= ln(medfaminc2000/(cpi2000/100))

* Farm production
gen lnvfprod30   = log(vfprod30/(cpi30/100))
gen lnvfprod40   = log(vfprod40/(cpi40/100))
gen lnvfprod60   = log(vfprod60/(cpi60/100))
gen lnvfprod2000 = log(vfprod2000/(cpi2000/100))

*Foreign born
gen fb0=fbwmtot00 + fbwftot00 
gen fb10=fbwtot10
gen fb20=fbwmtot20 + fbwftot20 
gen fb30=fbwmtot + fbwftot

gen fbshr20=fb20/(wmtot20 + wftot20)
gen fbshr30=fb30/(wmtot + wftot)

*Housing Values/Rents

ren medrnt30_NHGIS medrnt30
ren medhsval30_NHGIS medhsval30
ren var88_county72 medhsval70
ren var89_county72 medrnt70

foreach var in medhsval medrnt{
	foreach yr in 30 40 50 60 70 80 90 2000{
		cap replace `var'`yr'=`var'`yr'/(cpi`yr'/100)
	}
}

*various
drop other60

*drop counties experiencing big changes in area
gen d=(b1_lnd01_county00 - area)/(b1_lnd01_county00 + area)/2
drop if abs(d)>.03
replace area=(b1_lnd01_county00 + area)/2



******************************
* standardize variable names *
******************************


ren emp00 emp0
ren manuf_jobs_00 manuf_jobs_0

foreach yr in 0 10 20 30 40 50 60 70 80 90 2000{
	cap drop manuf`yr'
	cap drop agr`yr'
	ren manuf_jobs_`yr' manuf`yr'
	cap ren ag_jobs`yr' agr`yr'
}




foreach yr in 0 10 20 30 60 2000{
	gen other`yr'=emp`yr'-agr`yr'-manuf`yr'
}


*********************
*    Make Share     *
*********************

foreach var in manuf agr{
	foreach yr in 0 10 20 30 40 60 70 2000{
		cap gen `var'shr`yr'=`var'`yr'/emp`yr'
	}
}

********************
* prepare outcomes *
********************


foreach var in pop emp house wage twage agr manuf other medhsval medrnt fb {
	cap gen ln`var'0=ln(`var'0)
	cap gen ln`var'10=ln(`var'10)
	cap gen ln`var'20=ln(`var'20)
	cap gen ln`var'30=ln(`var'30)
	cap gen ln`var'40=ln(`var'40)
	cap gen ln`var'50=ln(`var'50)
	cap gen ln`var'60=ln(`var'60)
	cap gen ln`var'70=ln(`var'70)
	cap gen ln`var'80=ln(`var'80)
	cap gen ln`var'90=ln(`var'90)
	cap gen ln`var'2000=ln(`var'2000)
}

drop lntwage20  lnmedhsval20 lnmedrnt20 //stata confuses 20 and 2000


*******************
*Generate Controls*
*******************
g urbshare0=popurb0/pop0
g urbshare10=popurb10/pop10
g urbshare20=popurb20/pop20
g urbshare30=popurb30/pop20
gen popdens0=pop0/b1_lnd01_county00


*fix zeros in covariate quantities

foreach var in agr manuf other{
	foreach yr in 10 20 30{
		cap replace ln`var'`yr'=0 if `var'`yr'<=0
		gen no`var'`yr'dum=`var'`yr'<=0
	}
}

ren mfgcap00 mfgcap0
gen lnmfgcap0=ln(mfgcap0)
replace lnmfgcap0=0 if mfgcap0==0



*fix wages

foreach var in wage twage{
	foreach yr in 20 30{
		cap replace ln`var'`yr'=-1 if `var'`yr'==.
		cap gen no`var'`yr'dum=`var'`yr'==.
	}
}



*transformations
foreach var of varlist elevmax elevrang popdens0 tillit10 tillit1020 tillit1010 retsales radiorep totunemp area{
	gen ln`var'=ln(`var')
}


center tmean* lnelevmax

foreach var of varlist c_lnelevmax white0 white20 white30 c_tmean* lnmanuf0 lnmanuf10 lnmanuf20 lnradiorep lnemp20 lnemp30 lnwage0 lnwage20 lnwage30 lntwage30 lnpop0 lnpop20 lnpop30 lnfaval0 lnfaval20 lnfaval30 tmean*{
	gen `var'sq=`var'^2
	gen `var'cub=`var'^3
}

gen popdifsq=(lnpop30-lnpop20)^2
gen empdifsq=(lnemp30-lnemp20)^2
gen manufdifsq=(lnmanuf20-lnmanuf0)^2
gen urbsharedifsq=(urbshare20-urbshare0)^2
gen whitedifsq=(white20-white0)^2
gen agrdifsq=(lnagr20-lnagr0)^2
gen wagedifsq=(lnwage20-lnwage0)^2
gen favaldifsq=(lnfaval20-lnfaval0)^2

gen agrshr30sq=agrshr30^2
gen agrshr20sq=agrshr20^2

gen lnagr20sq=lnagr20^2
gen lnagr0sq=lnagr0^2
gen fbdifsq=(lnfb20-lnfb0)^2

gen wagedif2sq=(lnwage30-lnwage20)^2
gen tmean_jan_jul=tmean_jan*tmean_jul

gen pctil20=tillit1020/t10tot20
gen pctil30=tillit10/t10tot
replace PRADIO=PRADIO/100
gen urate30=totunemp/(totunemp+emp30)


global X "lnelevmax lnelevrang lnarea lnpop20 lnpop20sq lnpop30 lnpop30sq popdifsq agrshr20 agrshr20sq agrshr30 agrshr30sq manufshr20 manufshr30 nowage20 nowage30 lnwage20 lnwage30 notwage30 lntwage30 lnemp20 lnemp30 urbshare20 urbshare30 lnfaval20 lnfaval30 lnmedhsval30 lnmedrnt30 white20 white20sq white30 white30sq pctil20 pctil30 PRADIO urate30 fbshr20 fbshr30"

run "~/geo/x_ols_JP_v10.ado" //use Juan Pablo's version of Conley (1999)
run "~/geo/bo2_JP.ado" //use Juan Pablo's version of O-B
global tables="~/latex/tables_raw"



*******************************************
*******************************************
* Propensity score 
*******************************************
*******************************************

sum $X

logit tva $X, cluster(state)
predict phat

keep if e(sample)==1

sum phat if tva==1, det
sum phat if tva==0, det
local cut=r(p25)

tab tva



gen w=pop50 //weights


drop if phat<`cut'


*****************************************
*****************************************
* 	Reshape Data
*****************************************
*****************************************
tab region, gen(regdum)

foreach yr in 20 30 40 50 60 70 80 90 2000{
	cap gen agrdens`yr'=agr`yr'/b1_lnd01_county00
	cap gen lnagrdens`yr'=ln(agrdens`yr')
	gen manufdens`yr'=manuf`yr'/b1_lnd01_county00
	gen lnmanufdens`yr'=ln(manuf`yr'/b1_lnd01_county00)
}

*quantile knots for spline components
local knots "60 85" 
*local knots "33 66" //uncomment for alternate knots

centile manufdens80, c(`knots')
local cut1=r(c_1)
local cut2=r(c_2)

mkspline manufdens30_1 `cut1' manufdens30_2 `cut2' manufdens30_3 = manufdens30
mkspline manufdens40_1 `cut1' manufdens40_2 `cut2' manufdens40_3 = manufdens40
mkspline manufdens50_1 `cut1' manufdens50_2 `cut2' manufdens50_3 = manufdens50
mkspline manufdens60_1 `cut1' manufdens60_2 `cut2' manufdens60_3 = manufdens60
mkspline manufdens70_1 `cut1' manufdens70_2 `cut2' manufdens70_3 = manufdens70
mkspline manufdens80_1 `cut1' manufdens80_2 `cut2' manufdens80_3 = manufdens80
mkspline manufdens90_1 `cut1' manufdens90_2 `cut2' manufdens90_3 = manufdens90

centile lnmanufdens80, c(`knots')
local cut1=r(c_1)
local cut2=r(c_2)

mkspline lnmanufdens30_1 `cut1' lnmanufdens30_2 `cut2' lnmanufdens30_3 = lnmanufdens30
mkspline lnmanufdens40_1 `cut1' lnmanufdens40_2 `cut2' lnmanufdens40_3 = lnmanufdens40
mkspline lnmanufdens50_1 `cut1' lnmanufdens50_2 `cut2' lnmanufdens50_3 = lnmanufdens50
mkspline lnmanufdens60_1 `cut1' lnmanufdens60_2 `cut2' lnmanufdens60_3 = lnmanufdens60
mkspline lnmanufdens70_1 `cut1' lnmanufdens70_2 `cut2' lnmanufdens70_3 = lnmanufdens70
mkspline lnmanufdens80_1 `cut1' lnmanufdens80_2 `cut2' lnmanufdens80_3 = lnmanufdens80
mkspline lnmanufdens90_1 `cut1' lnmanufdens90_2 `cut2' lnmanufdens90_3 = lnmanufdens90


*for graphing
local thresh1 `cut1'
local thresh2 `cut2'

centile lnmanufdens80, c(5 10 15 20 25 30 35 40 45 50)
local cutg5=r(c_1)
local cutg10=r(c_2)
local cutg15=r(c_3)
local cutg20=r(c_4)
local cutg25=r(c_5)
local cutg30=r(c_6)
local cutg35=r(c_7)
local cutg40=r(c_8)
local cutg45=r(c_9)
local cutg50=r(c_10)

centile lnmanufdens80, c(55 60 65 70 75 80 85 90 95)
local cutg55=r(c_1)
local cutg60=r(c_2)
local cutg65=r(c_3)
local cutg70=r(c_4)
local cutg75=r(c_5)
local cutg80=r(c_6)
local cutg85=r(c_7)
local cutg90=r(c_8)
local cutg95=r(c_9)

*make spline terms
forvalues i=1/3{
	gen lDmanufdens_`i'=manufdens90_`i'-manufdens80_`i'
	gen l2Dmanufdens_`i'=manufdens80_`i'-manufdens70_`i'
	gen l3Dmanufdens_`i'=manufdens70_`i'-manufdens60_`i'
	gen l4Dmanufdens_`i'=manufdens60_`i'-manufdens50_`i'
	gen l5Dmanufdens_`i'=manufdens50_`i'-manufdens40_`i'
	gen l6Dmanufdens_`i'=manufdens40_`i'-manufdens30_`i'

	gen lDlnmanufdens_`i'=lnmanufdens90_`i'-lnmanufdens80_`i'
	gen l2Dlnmanufdens_`i'=lnmanufdens80_`i'-lnmanufdens70_`i'
	gen l3Dlnmanufdens_`i'=lnmanufdens70_`i'-lnmanufdens60_`i'
	gen l4Dlnmanufdens_`i'=lnmanufdens60_`i'-lnmanufdens50_`i'
	gen l5Dlnmanufdens_`i'=lnmanufdens50_`i'-lnmanufdens40_`i'
	gen l6Dlnmanufdens_`i'=lnmanufdens40_`i'-lnmanufdens30_`i'
}

gen lDmanufdens=manufdens90-manufdens80
gen l2Dmanufdens=manufdens80-manufdens70
gen l3Dmanufdens=manufdens70-manufdens60
gen l4Dmanufdens=manufdens60-manufdens50
gen l5Dmanufdens=manufdens50-manufdens40
gen l6Dmanufdens=manufdens40-manufdens30

gen lDlnmanufdens=lnmanufdens90-lnmanufdens80
gen l2Dlnmanufdens=lnmanufdens80-lnmanufdens70
gen l3Dlnmanufdens=lnmanufdens70-lnmanufdens60
gen l4Dlnmanufdens=lnmanufdens60-lnmanufdens50
gen l5Dlnmanufdens=lnmanufdens50-lnmanufdens40
gen l6Dlnmanufdens=lnmanufdens40-lnmanufdens30




gen D00=lnmanufdens2000-lnmanufdens90
gen D90=lnmanufdens90-lnmanufdens80
gen D80=lnmanufdens80-lnmanufdens70
gen D70=lnmanufdens70-lnmanufdens60
gen D60=lnmanufdens60-lnmanufdens50
gen D50=lnmanufdens50-lnmanufdens40


gen W00=lnwage2000-lnwage90
gen W90=lnwage90-lnwage80
gen W80=lnwage80-lnwage70
gen W70=lnwage70-lnwage60
gen W60=lnwage60-lnwage50
gen W50=lnwage50-lnwage40


corr D00 D90 D80 D70 D60 D50


**spline the instrument**
gen dens82=emp82/b1_lnd01_county00
gen dens77=emp77/b1_lnd01_county00
gen dens72=emp72/b1_lnd01_county00
gen dens67=emp67/b1_lnd01_county00
gen dens63=emp63/b1_lnd01_county00
gen dens54=emp54/b1_lnd01_county00
gen dens58=emp58/b1_lnd01_county00
gen dens47=emp47/b1_lnd01_county00


centile dens82, c(`knots')
local cut1=r(c_1)
local cut2=r(c_2)

mkspline dens82_1 `cut1' dens82_2 `cut2' dens82_3 = dens82
mkspline dens77_1 `cut1' dens77_2 `cut2' dens77_3 = dens77
mkspline dens72_1 `cut1' dens72_2 `cut2' dens72_3 = dens72
mkspline dens67_1 `cut1' dens67_2 `cut2' dens67_3 = dens67
mkspline dens63_1 `cut1' dens63_2 `cut2' dens63_3 = dens63
mkspline dens54_1 `cut1' dens54_2 `cut2' dens54_3 = dens54
mkspline dens58_1 `cut1' dens58_2 `cut2' dens58_3 = dens58
mkspline dens47_1 `cut1' dens47_2 `cut2' dens47_3 = dens47

forvalues i=1/3{
  gen D82_`i'= dens82_`i'-dens72_`i'
  gen D77_`i'= dens77_`i'-dens67_`i'
  gen D72_`i'= dens72_`i'-dens63_`i'
  gen D67_`i'= dens67_`i'-dens58_`i'
  gen D63_`i'= dens63_`i'-dens54_`i'
  gen D54_`i'= dens54_`i'-dens47_`i'
}


gen lndens82=ln(dens82)
gen lndens77=ln(dens77)
gen lndens72=ln(dens72)
gen lndens67=ln(dens67)
gen lndens63=ln(dens63)
gen lndens54=ln(dens54)
gen lndens58=ln(dens58)
gen lndens47=ln(dens47)


centile lndens82, c(`knots')
local cut1=r(c_1)
local cut2=r(c_2)

mkspline lndens82_1 `cut1' lndens82_2 `cut2' lndens82_3 = lndens82
mkspline lndens77_1 `cut1' lndens77_2 `cut2' lndens77_3 = lndens77
mkspline lndens72_1 `cut1' lndens72_2 `cut2' lndens72_3 = lndens72
mkspline lndens67_1 `cut1' lndens67_2 `cut2' lndens67_3 = lndens67
mkspline lndens63_1 `cut1' lndens63_2 `cut2' lndens63_3 = lndens63
mkspline lndens54_1 `cut1' lndens54_2 `cut2' lndens54_3 = lndens54
mkspline lndens58_1 `cut1' lndens58_2 `cut2' lndens58_3 = lndens58
mkspline lndens47_1 `cut1' lndens47_2 `cut2' lndens47_3 = lndens47

forvalues i=1/3{
  gen LD82_`i'= lndens82_`i'-lndens72_`i'
  gen LD77_`i'= lndens77_`i'-lndens67_`i'
  gen LD72_`i'= lndens72_`i'-lndens63_`i'
  gen LD67_`i'= lndens67_`i'-lndens58_`i'
  gen LD63_`i'= lndens63_`i'-lndens54_`i'
  gen LD54_`i'= lndens54_`i'-lndens47_`i'
}



**Detour: Get Deltas Before Completing Reshape**

**Preferred
*these coefficients come from Table 5
local b1=0.097
local b2=0.042
local b3=0.001


gen R4060=lnmanuf60-lnmanuf40-`b1'*(manufdens50_1-manufdens30_1)-`b2'*(manufdens50_2-manufdens30_2)-`b3'*(manufdens50_3-manufdens30_3) + 1.5*(lnwage60-lnwage40)
gen R6080=lnmanuf80-lnmanuf60-`b1'*(manufdens70_1-manufdens50_1)-`b2'*(manufdens70_2-manufdens50_2)-`b3'*(manufdens70_3-manufdens50_3) + 1.5*(lnwage80-lnwage60)
gen R802000=lnmanuf2000-lnmanuf80-`b1'*(manufdens90_1-manufdens70_1)-`b2'*(manufdens90_2-manufdens70_2)-`b3'*(manufdens90_3-manufdens70_3) + 1.5*(lnwage2000-lnwage80)


ivreg2 R4060 tva $X [aw=pop50], cluster(state) partial($X)
eststo delt50
ivreg2 R6080 tva $X  [aw=pop50], cluster(state) partial($X)
eststo delt70
ivreg2 R802000 tva $X  [aw=pop50], cluster(state) partial($X)
eststo delt90


*These come from Table 6

local b1=.4430518
local b2=.4562646
local b3=.4655453
 
gen LR4060=lnmanuf60-lnmanuf40-`b1'*(lnmanufdens50_1-lnmanufdens30_1)-`b2'*(lnmanufdens50_2-lnmanufdens30_2)-`b3'*(lnmanufdens50_3-lnmanufdens30_3) + 1.5*(lnwage60-lnwage40)
gen LR6080=lnmanuf80-lnmanuf60-`b1'*(lnmanufdens70_1-lnmanufdens50_1)-`b2'*(lnmanufdens70_2-lnmanufdens50_2)-`b3'*(lnmanufdens70_3-lnmanufdens50_3) + 1.5*(lnwage80-lnwage60)
gen LR802000=lnmanuf2000-lnmanuf80-`b1'*(lnmanufdens90_1-lnmanufdens70_1)-`b2'*(lnmanufdens90_2-lnmanufdens70_2)-`b3'*(lnmanufdens90_3-lnmanufdens70_3) + 1.5*(lnwage2000-lnwage80)

ivreg2 LR4060 tva $X  [aw=pop50], cluster(state) partial($X)
eststo delt50L
ivreg2 LR6080 tva $X   [aw=pop50], cluster(state) partial($X)
eststo delt60L
ivreg2 LR802000 tva $X  [aw=pop50], cluster(state) partial($X)
eststo delt80L

esttab delt* using $tables/delts.csv, replace keep(tva) b(3) se(3) star(* 0.1 ** 0.05 *** 0.01) mlab(40-60 60-80 80-2000 "40-60 logs" "60-80 logs" "80-2000 logs")


*convert to long form (Pool changes)
ren D00 D100
ren W00 W100

ren D82_1 Z1_100
ren D82_2 Z2_100
ren D82_3 Z3_100

ren D72_1 Z1_90
ren D72_2 Z2_90
ren D72_3 Z3_90

ren D63_1 Z1_80
ren D63_2 Z2_80
ren D63_3 Z3_80

ren D54_1 Z1_70
ren D54_2 Z2_70
ren D54_3 Z3_70


ren LD82_1 LZ1_100
ren LD82_2 LZ2_100
ren LD82_3 LZ3_100

ren LD72_1 LZ1_90
ren LD72_2 LZ2_90
ren LD72_3 LZ3_90

ren LD63_1 LZ1_80
ren LD63_2 LZ2_80
ren LD63_3 LZ3_80

ren LD54_1 LZ1_70
ren LD54_2 LZ2_70
ren LD54_3 LZ3_70


ren lDmanufdens_1 X1_100
ren lDmanufdens_2 X2_100
ren lDmanufdens_3 X3_100

ren l2Dmanufdens_1 X1_90
ren l2Dmanufdens_2 X2_90
ren l2Dmanufdens_3 X3_90

ren l3Dmanufdens_1 X1_80
ren l3Dmanufdens_2 X2_80
ren l3Dmanufdens_3 X3_80

ren l4Dmanufdens_1 X1_70
ren l4Dmanufdens_2 X2_70
ren l4Dmanufdens_3 X3_70

ren l5Dmanufdens_1 X1_60
ren l5Dmanufdens_2 X2_60
ren l5Dmanufdens_3 X3_60

ren l6Dmanufdens_1 X1_50
ren l6Dmanufdens_2 X2_50
ren l6Dmanufdens_3 X3_50


ren lDlnmanufdens_1 LX1_100
ren lDlnmanufdens_2 LX2_100
ren lDlnmanufdens_3 LX3_100

ren l2Dlnmanufdens_1 LX1_90
ren l2Dlnmanufdens_2 LX2_90
ren l2Dlnmanufdens_3 LX3_90

ren l3Dlnmanufdens_1 LX1_80
ren l3Dlnmanufdens_2 LX2_80
ren l3Dlnmanufdens_3 LX3_80

ren l4Dlnmanufdens_1 LX1_70
ren l4Dlnmanufdens_2 LX2_70
ren l4Dlnmanufdens_3 LX3_70

ren l5Dlnmanufdens_1 LX1_60
ren l5Dlnmanufdens_2 LX2_60
ren l5Dlnmanufdens_3 LX3_60

ren l6Dlnmanufdens_1 LX1_50
ren l6Dlnmanufdens_2 LX2_50
ren l6Dlnmanufdens_3 LX3_50


ren lDmanufdens X_100
ren l2Dmanufdens X_90
ren l3Dmanufdens X_80
ren l4Dmanufdens X_70
ren l5Dmanufdens X_60
ren l6Dmanufdens X_50

ren lDlnmanufdens LX_100
ren l2Dlnmanufdens LX_90
ren l3Dlnmanufdens LX_80
ren l4Dlnmanufdens LX_70
ren l5Dlnmanufdens LX_60
ren l6Dlnmanufdens LX_50

*ren lnmanufdens40 lnmanufdens_40
ren lnmanufdens50 lnmanufdens_50
ren lnmanufdens60 lnmanufdens_60
ren lnmanufdens70 lnmanufdens_70
ren lnmanufdens80 lnmanufdens_80
ren lnmanufdens90 lnmanufdens_90
ren lnmanufdens2000 lnmanufdens_100

keep *Z* W* D100 D90 D80 D70 D60 D50 *X1* *X2* *X3* *X_* fips state region $X lnmanufdens40 pop80 pop50 b1_lnd01_county00 tva lnmanufdens_* L* b1_lnd01_county00
drop *67* *77*
reshape long W Z1_ Z2_ Z3_ LZ1_ LZ2_ LZ3_ D X1_ X2_ X3_ LX1_ LX2_ LX3_ X_ lnmanufdens_ LX, i(fips state region lnmanufdens40 pop80 pop50 b1_lnd01_county00 tva b1_lnd01_county00) j(t)
sum



tsset fips t, delta(10)

gen yr90=t==90
gen yr80=t==80
gen yr70=t==70


tab region, gen(regdum)


**Calibrated lhs variable
gen R=D+1.5*W


***ESTIMATION***

**spline in levels**
ivreg2 R X1 X2 X3 yr70 yr80 yr90 $X tva* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90) ffirst
eststo spline1
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)

ivreg2 R X1 X2 X3 yr70 yr80 yr90 $X tva* lnmanufdens40* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90 lnmanufdens40*) ffirst
eststo spline2
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)

ivreg2 R X1 X2 X3 yr70 yr80 yr90 $X tva* lnmanufdens40* regdum* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90 lnmanufdens40* regdum*) ffirst
eststo spline3
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)


ivreg2 R (X1 X2 X3 = Z1 Z2 Z3 ) yr70 yr80 yr90 $X tva* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90) ffirst
eststo spline4
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)

ivreg2 R (X1 X2 X3= Z1 Z2 Z3) yr70 yr80 yr90 $X tva* lnmanufdens40* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90 lnmanufdens40*) ffirst
eststo spline5
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)

ivreg2 R (X1 X2 X3= Z1 Z2 Z3) yr70 yr80 yr90 $X tva* lnmanufdens40* regdum* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90 lnmanufdens40* regdum*) ffirst
eststo spline6
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)


esttab spline1 spline2 spline3 spline4 spline5 spline6 using $tables/g_terciles_opt.csv, replace cell( b(fmt(a2)) se(par fmt(a2)) first[APF](par([ ]) fmt(2) ) ) stat(equal all N, fmt(%6.0g))

**now in logs**
ivreg2 R LX1 LX2 LX3 yr70 yr80 yr90 $X tva* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90) ffirst
eststo spline1
test LX1=LX2=LX3
estadd scalar equal=r(p)
test LX1 LX2 LX3
estadd scalar all=r(p)

ivreg2 R LX1 LX2 LX3 yr70 yr80 yr90 $X tva* lnmanufdens40* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90 lnmanufdens40*) ffirst
eststo spline2
test LX1=LX2=LX3
estadd scalar equal=r(p)
test LX1 LX2 LX3
estadd scalar all=r(p)

ivreg2 R LX1 LX2 LX3 yr70 yr80 yr90 $X tva* lnmanufdens40* regdum* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90 lnmanufdens40* regdum*) ffirst
eststo spline3
test LX1=LX2=LX3
estadd scalar equal=r(p)
test LX1 LX2 LX3
estadd scalar all=r(p)


ivreg2 R (LX1 LX2 LX3 = LZ1 LZ2 LZ3 ) yr70 yr80 yr90 $X tva* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90) ffirst
eststo spline4
test LX1=LX2=LX3
estadd scalar equal=r(p)
test LX1 LX2 LX3
estadd scalar all=r(p)

ivreg2 R (LX1 LX2 LX3= LZ1 LZ2 LZ3) yr70 yr80 yr90 $X tva* lnmanufdens40* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90 lnmanufdens40*) ffirst
eststo spline5
test LX1=LX2=LX3
estadd scalar equal=r(p)
test LX1 LX2 LX3
estadd scalar all=r(p)

ivreg2 R (LX1 LX2 LX3= LZ1 LZ2 LZ3) yr70 yr80 yr90 $X tva* lnmanufdens40* regdum* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90 lnmanufdens40* regdum*) ffirst
eststo spline6
test LX1=LX2=LX3
estadd scalar equal=r(p)
test LX1 LX2 LX3
estadd scalar all=r(p)


esttab spline1 spline2 spline3 spline4 spline5 spline6 using $tables/g_terciles_opt_logs.csv, replace cell( b(fmt(3)) se(par fmt(3)) first[APF](par([ ]) fmt(2) ) ) stat(equal all N, fmt(%6.0g))


**sensitivity**
gen R1=D+W
gen R2=D+2*W

ivreg2 R1 (X1 X2 X3 = Z1 Z2 Z3 ) yr70 yr80 yr90 $X tva* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90) ffirst
eststo spline1
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)

ivreg2 R1 (X1 X2 X3 = Z1 Z2 Z3 ) yr70 yr80 yr90 $X tva* lnmanufdens40* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90 lnmanufdens40*) ffirst
eststo spline2
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)

ivreg2 R1 (X1 X2 X3 = Z1 Z2 Z3 ) yr70 yr80 yr90 $X tva* lnmanufdens40* regdum* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90 lnmanufdens40* regdum*) ffirst
eststo spline3
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)


ivreg2 R2 (X1 X2 X3 = Z1 Z2 Z3 ) yr70 yr80 yr90 $X tva* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90) ffirst
eststo spline4
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)

ivreg2 R2 (X1 X2 X3= Z1 Z2 Z3) yr70 yr80 yr90 $X tva* lnmanufdens40* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90 lnmanufdens40*) ffirst
eststo spline5
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)

ivreg2 R2 (X1 X2 X3= Z1 Z2 Z3) yr70 yr80 yr90 $X tva* lnmanufdens40* regdum* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90 lnmanufdens40* regdum*) ffirst
eststo spline6
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)


esttab spline1 spline2 spline3 spline4 spline5 spline6 using $tables/g_terciles_opt_sens.csv, replace cell( b(fmt(a2)) se(par fmt(a2)) first[APF](par([ ]) fmt(2) ) ) stat(equal all N, fmt(%6.0g))



**estimate wage coeff**
ivreg2 D W X1 X2 X3 yr70 yr80 yr90 $X tva* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90) ffirst
eststo spline1
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)

ivreg2 D W X1 X2 X3 yr70 yr80 yr90 $X tva* lnmanufdens40* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90 lnmanufdens40*) ffirst
eststo spline2
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)

ivreg2 D W X1 X2 X3 yr70 yr80 yr90 $X tva* lnmanufdens40* regdum* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90 lnmanufdens40* regdum*) ffirst
eststo spline3
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)


ivreg2 D W (X1 X2 X3 = Z1 Z2 Z3 ) yr70 yr80 yr90 $X tva* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90) ffirst
eststo spline4
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)

ivreg2 D W (X1 X2 X3= Z1 Z2 Z3) yr70 yr80 yr90 $X tva* lnmanufdens40* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90 lnmanufdens40*) ffirst
eststo spline5
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)

ivreg2 D W (X1 X2 X3= Z1 Z2 Z3) yr70 yr80 yr90 $X tva* lnmanufdens40* regdum* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90 lnmanufdens40* regdum*) ffirst
eststo spline6
test X1=X2=X3
estadd scalar equal=r(p)
test X1 X2 X3
estadd scalar all=r(p)


esttab spline1 spline2 spline3 spline4 spline5 spline6 using $tables/g_terciles_opt_freewage.csv, replace cell( b(fmt(a2)) se(par fmt(a2)) first[APF](par([ ]) fmt(2) ) ) stat(equal all N, fmt(%6.0g))




**info for graphs**
qui ivreg2 R (LX1 LX2 LX3= LZ1 LZ2 LZ3) yr70 yr80 yr90 $X tva* [w=pop50] if t>=70, cluster(state) partial($X yr70 yr80 yr90) ffirst

di `thresh1' `thresh2'
mat A=J(19,3,.)
local j=1
forvalues i=5(5)95{
  local newcut `cutg`i''
  di `newcut'
  
  if `newcut' <= `thresh1' {
    lincom LX1*`newcut'
  }
  if `newcut' > `thresh1' {
    lincom LX1*`thresh1' + LX2*(`newcut'-`thresh1')
  }
  if `newcut' > `thresh2' {
    lincom LX1*`thresh1' + LX2*(`thresh2'-`thresh1') + LX3*(`newcut'-`thresh2')
  }

  mat A[`j',1]=`newcut'
  mat A[`j',2]=r(estimate)
  mat A[`j',3]=r(se)
  local ++j
  }
  
  svmat A
  gen ub=A2+1.96*A3
  gen lb=A2-1.96*A3
  list A2 ub lb in 1/20

ren A1 x_est
ren A2 g_est
ren ub ub_est
ren lb lb_est
drop A*

keep x_est g_est ub_est lb_est
save gfuncgraph, replace



